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O ' Abstract 
(N 

^ , The paper focuses on minimum mean square error (MMSE) Bayesian estimation for a Gaussian 

■ source impaired by additive Middleton's Class-A impulsive noise. In order to reduce the implementation 

' complexity associated with the expression of the optimum Bayesian estimator, the paper considers also 

(N ■ 

two popular suboptimal estimators, which are the soft-limiter and the blanker. The optimum Bayesian 



thresholds for such suboptimal estimators are obtained by iteratively solving fixed point equations. 
Connections with the maximum SNR estimators are also estabhshed. Theoretical expressions for the MSE 
and the SNR of the suboptimal estimators are also derived, and the results confirmed by simulations. 
Noteworthy, the results can be applied to any noise characterized by a Gaussian mixture distribution. 

I. INTRODUCTION 

Interference and noise with impulsive non-Gaussian distributions may impair the performance 
of several systems including communications, controls, sensors and so forth. Middleton has proposed 
widely accepted canonical models for interference lll-lJl, which are capable to characterize "intelligent" 
(e.g., information bearing), as well as "non-intelligent" (e.g., natural or man-made) noises. These models 
^ • distinguish impulsive interference confined within the observer bandwidth or, equivalently, which does 

not introduce transitory (i.e., Class-A noise), from wideband interference inducing a transitory in the 
observing system (i.e., Class-B and Class-C noises) ||4l. Although Middleton's noise models were widely 
investigated to identify the interference behavior |[Tl-|[5l, to estimate their canonical parameters ll6l- |[T0l . 
and to detect finite alphabets in digital communications |[m - |[T6l . to the best of the author knowledge, 
results for the optimum Bayesian estimator (OBE) of Gaussian sources in Class-A impulsive noise are 
still lacking or, at least, not widely acknowledged. Actually, such OBE estimator would be of valuable 
help in those applications where the quantity of interest can be modeled, or approximated, by a Gaussian 
probability density function (pdf). This happens, thanks to the central limit theorem (CLT) ifTTl . when the 
quantity of interest is generated by the superposition of several other non-dominating random quantities. 
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This is the case, for instance, of a communication system when several signals, possibly belonging to 
different users, are superimposed in time and multiplexed in another domain, such as the frequency or the 
code domain. Well know examples are multicarrier systems, such as asymmetric digital subscriber lines 
(ADSL) ifTSl and orthogonal frequency division multiplexing (OFDM) ||T9l . or code division multiple 
access (CDMA) systems |[T9l . More generally, any multiple-input single-output (MISO) system with a 
high number of inputs, can be approximately modeled this way. 

The first aim of the paper is to derive the minimum mean squared error (MMSE) OBE for a scalar 
Gaussian source impaired by Middleton's Class-A canonical noise. Noteworthy, such an OBE could be 
useful also as a preprocessing stage for those estimation and detection algorithms that are designed 
under AWGN hypotheses and, as such, are not robust to impulsive noises lITTl . 1201 . Although the paper 
derives the analytical expression of the OBE in a closed form, its use may be restricted in some practical 
applications due to complexity constraints. In these cases, protection from impulsive noise may be granted 
by simpler suboptimal devices that are robust to high noise peaks. A possibility is to resort either to a 
blanking-nonlinearity (BN) that nulls out the received signal when it overpasses a given threshold, or, 
alternatively, to a soft-limiter (SL) that simply clips the signal when it overpasses the threshold. In both 
cases, only the blanking or the clipping thresholds have to be optimized in a MMSE Bayesian sense. This 
way, although both the Bayesian BN estimator (BNE) or the Bayesian SL estimator (SLE) are suboptimal 
with respect to the OBE, they greatly reduce the implementation complexity, because they only request 
to store a scalar value (e.g., the threshold), rather than the OBE curve, for the set of source and noise 
statistical parameters of interest. 

Actually, although the SLE and BNE are suboptimal with respect to the OBE, the derivation of the 
optimum Bayesian thresholds is analytically much harder than the computation of the OBE expression. 
Specifically, the computation of the optimum thresholds can be formulated in both cases as the solution 
of a fixed-point problem 11211 . which is proved to always admit a solution, obtainable by standard iterative 
approaches with fast convergence. 

The comparison of the shape of the OBE curve with the simplified SLE and BNE, can intuitively 
illuminate whether the BNE or the SLE is the best simplified strategy. The MSE associated to the SLE 
and the BLE, is the quantitative parameter to choose among the two, and it highlights that the best choice 
strictly depends on the statistical characteristics of the impulsive noise, which are summarized by the 
average signal-to-noise power ratio (SNR), the peakness, the average number of emitting noise sources, 
etc.. 

A different criterion to set the optimal threshold of a BN, has been proposed in EOl . which is based 



3 



on the maximization of the SNR at the estimator output. The author in ll20l exploits an extension of 
the Bussgang theorem |[22l . ifTTl to complex Gaussian mixtures at the input of a nonlinearity (i.e., of 
the BN). Similar results has been derived also for a SL in a two-component Gaussian mixture noise in 
||23]| . Note that, while MMSE and maximum-SNR (MSNR) are equivalent in pure AWGN scenarios |l24l. 
where the MMSE estimator is linear, this is not the case when the noise is a Gaussian-mixture, which 
leads to a non-linear MMSE estimator. This paper establishes also the connection between the MMSE 
and the MSNR criteria, at the output of the non-linear estimator. Note that, while ll20l and |[23l consider 
a complex Gaussian source, and derive the optimum MSNR thresholds for BN and SL that operate on 
the complex signal envelope, this paper concentrates on real Gaussian sources impaired by Middleton's 
Class-A noises. Thus, the optimal MSNR thresholds, derived in this paper for the BN and the SL, are 
different from those derived in |[20l and |[23l . If, anyway, the source is a complex Gaussian random 
variable, the results in this paper can be applied to SL and BN that separately operate on each Cartesian 
components of the complex source. 

Due to the fact that the MMSE and the MSNR approaches are not equivalent, they lead to suboptimal 
estimators with different thresholds: this paper shows when the two thresholds are similar, and when, 
conversely, they are different. Whether it is better to maximize the SNR or minimize the MSE, depends 
on the specific appUcation and design constraint. This is not however the subject of the paper. The paper is 
organized as follows: section HIl introduces the system model, while the OBE, SLE, and BNE are derived 
in sections |lIIJ |IVl and |Vl respectively. Successively, section |Vl] establishes the connection between the 
MMSE and MSNR criteria, and derives the MSE expressions for the SLE and the BNE. Finally section 
IVIII is dedicated to computer simulations that confirm the theoretical findings, while the conclusions are 
drawn in the last section. 



Let's consider a real source x with average power a\, and zero-mean Gaussian pdf, as expressed by 



which is impaired by a Class-A impulsive noise n, with average power cj^, as shown in Fig. [T] and 
summarized by 



IL System Model 



fx{x) = G{x;a\) 



1 
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(1) 



y = X + n. 



(2) 
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Fig. 1 : System model 



The Class- A impulsive noise subsumes also the presence of a background zero-mean thermal AWGN nj, 
with average power erf . Specifically, the impulsive noise n has a pdf expressed by 

oo oo QO n 2 

fN{n) = fimfmin) = J2 PraG{n-al) = J2 Tf^e""^' (3) 

m=0 m=0 m=0 v m 

which is clearly a Gaussian-mixture, where the weights (3m = e~^A"^ /ml represent the Poisson- 

distributed probability that m noise sources simultaneously contribute to the impulsive event 111, 

The power associated with the noise emission from the m noise sources, is expressed by 

2 m/A + T 2 cr? 2 

= [^r^ (^N = + f^t , (4) 

where o"^ = E{n'^} = aj + a^, aj represents the impulsive part of the noise power, T = cr^/aj is 

oo 

the power ratio among the AWGN and the impulsive part of the noise n, and A = E{m} = J2 i^T-Pm 

m=0 

represents the average number of impulsive sources that are simultaneously active. Note that E{-} is used 
throughout the paper for statistical expectation, and Ex{-} is used to make explicit that the expectation 
is computed with respect to the pdf of the random variable X. 

The overall Class-A model is typically identified by the three canonical parameters. A, T and, cr^, 
which summarize the statistical structure of the noise. Specifically, low values of A identify rare and 
highly peaked impulsive sources, and, conversely, high values of A, by a CLT argument, makes the 
impulsive noise more similar to an AWGN. The interested readers are redirected to 13, H, 0, and 
references therein, for further details on the Class-A model, and to ll6l- |[T0l . for the estimation of the 
canonical parameters A, T, and aj^. This paper assumes that the canonical parameters have been perfectly 
estimated by the observing system. 
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III. Optimum Bayesian Estimator (OBE) 
The optimum non-linear MMSE Bayesian estimator of x given the observed data y, is expressed by 

m 

+00 

xoB¥.{y) = Ex\y{x\v} = j xfx\Y{x;y)dx, (5) 

—00 

where fx\Y{x; y) represents the conditional pdf of the source x for a given observation y. By exploiting 

Bayes rules for conditional pdf?,, ^ can be expressed as 

00 

ioBE (2/) = -r-TT / ^fY\x {y\x) fx{x)dx. (6) 

—00 

Due to the fact that the impulsive noise n is independent of x, it is well known lITTl that the probability 
density function of the observed value y in Q is given by the linear convolution of the Gaussian pdf 

fx{y) in O with the Middleton-A pdf fN{y) in ([3]), as expressed by 

00 

fy (y) = E l^rnGiy; al) * Gx{y; (7) 

m=0 

where * stands for the linear convolution operator. 

Thus, the convolution of two Gaussian pdf?, generates a new Gaussian pdf fTTI], with variance equal 
to the sum of the two original variances, and, consequently, equation ([7]) becomes 



00 00 n 

fy (y) = E PmG{y; al + ai) = ^ , e 



m=0 m=0 y 27r(cj2^ + (T 

Observing (|2]|, it is also evident that fY\x{y] x) in Q, is expressed by 



4. 



(8) 



fY\x{y:x) = fN{y-x)=Y.^=e (9) 

Plugging © and ([TJ in (|6]l, the overall integral becomes 

00 00 

^OBE(y) = -L Xfx{x)fm {y - X) dx 

_ 1 

which contains the convolution between the m-th Gaussian pdf fm{x) and p{x) = xfx{x). The con- 
volution can be computed in the Fourier domain by exploiting the properties of the Fourier transform 
operator, as detailed in Appendix A, leading to 



m=0 —00 

^ l^mlpiy)* fm{y)], 

m=0 



(10) 



xoBE{y) = — y- (11) 

m=0 ^ "> ' X/ 
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Equation ([TT]) highlights how the OBE depends on the source average power aj^, and the noise canonical 
parameters A, T, and aj^, through (3m and in ©. The input-output characteristic of the OBE is 
plotted in Fig. 2(a)| and Fig. |2(b)| for several values of the parameter A that controls the peakness of 
the impulsive noise fSl; it is evident that for very high values of A, when /Ar(n) tends to a zero-mean 
Gaussian pdf, the OBE tends to the well known linear-MMSE estimator |25|, expressed by 

-(/in). ^ E{xy} a\ 

-'o.i{y) = ^^y = ^^^v- (12) 

Conversely, for lower values of A, when the noise is characterizes by rare and highly peaked impulses, 
the optimum estimator shows a highly non-linear nature, by roughly limiting (A = 0.1) or blanking 
(A = 0.001) the observed values y that overpass certain thresholds. Note that, the OBE expression in 
([TT]) . immediately applies also to any other Gaussian-mixture noise n, possibly different from the Class-A 

Af-l 

noise, with a finite number M of mixtures, such that /3m > 0, and Yl, f^m = 1- 

m=0 

IV. Bayesian Soft Limiter Estimator (SLE) 

The OBE rather involved analytical expression ([TTb could prevent its use in real-time applications that 
are subject to either memory or computational complexity constraints. This is especially true when the 
OBE analytical expression is requested to be adaptive with respect to changes of either the source or the 
noise statistical parameters (e.g., the average powers aj^ and cr^, and the noise peakness factors A, and 
T). For this reason, this section investigates a simple suboptimum estimator, namely the SL estimator 
shown in Fig. [3l which adds robustness to a system impaired by impulsive noise events, by clipping 
those values that exceed a given threshold a. Thus, the only parameter to optimize in the Bayesian sense 
is the clipping threshold a, which obviously would depend on the noise parameters A, T and, afj, as 
well as on the source power a^- Meaningfulness of such an optimization, which leads to the SLE is 
also suggested by the OBE shapes in Fig. |2(a)| and Fig. 2(b) which for certain noise parameters (e.g., 
^ = 0.1) resembles the SL estimator of Fig. [3l 

Note that, the output of the SL in Fig. [3] is expressed by a non-linear input-output characteristic 
S:sL = 9sL(y;a)- Thus, the SL estimation error csl depends on the selected threshold a, as well as on 
the statistical properties of the source x and the noise n. This is expressed by 

esL = hshix, n;a) = X- qsl^x + n; a), (13) 

where hsL{x,n;a) highlights that the estimation error csl is a non linear transformation of two random 
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Fig. 2: OBE for several values of A 



variables, x and n, as summarized by 

hsLix,n;a) 



X + a ,x + n< —a 
—n ,\x + n\ < a ■ 
X — a , x + n > a 



(14) 



The SLE estimator is defined by selecting ag™*^-* according to the MMSE criterion, as it follows 



a 



(mse) 
SL 



arg mm 

Q!e7^+ 



E 



{4.} 



arg mm 



E 



{hl^{x,n;a)^ 



(15) 
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Fig. 3: Soft Limiter estimator 



Thus, in order to find a^^ it is necessary to solve 

i^E{hl^{x,n;a)} = E U [hl{x,n; a)]} 

f (la) 1 

= E <2hsLix,n; a)h^^^ ' {x,n;a)> = 0, 

where h^^^"\x,n; a) stands for the first partial derivative of hsi^{x,n; a) with respect to a, and it leads 
to the following integral equation 

+ 00 +00 

hsL{x,n;a)h^^{°'\x,n;a)fx{x)fN{n)dxdn = 0. (17) 



By substituting ([T), Q, (fT4l) . and its partial derivative in ([TT] ). after tedious derivations detailed in 
Appendix B, it is proved that a^^ is the solution of the following fixed-point equation 

a = ^SL («) = 2cTx 33 ; , (18) 

1- E /3rrxerf(a/J2(a|+a2 



Appendix B proves that the solution of the fixed point equation (118]) always exists. Moreover, by 
conjecturing and proving by simulations that F^'^^'^\a) is a contraction mapping |[2TI before its first 
maximum, it is inferred that the solution is also unique. Thus, the succession a„+i = F^^'^^\an) 
converges to the exact fixed point solution a^™'^^^ when n goes to infinity lllll . Consequently, Og'^'^'^^can 
be numerically approximated by the following iterative algorithm 



Al: Iterative algorithm for optimal SL threshold 
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set ao = ^SL (0) ^nd n = 0; 

while \F^'^'^^\an) — a„| > e and n < um, 

j-,(mse) , X 
an+1 = -FsL ' 

n = n + 1; 

end 

(mse) 

set = a„. 



In algorithm Al, e represents the accuracy that is requested for the approximated solution to stop 
within a maximum number nMax of iterations. Obviously, other iterative numerical approaches can also 
be used to solve dTSl ). such as the Newton-Rapson method |[211 to find the root of F^^'^^'^ (a) — a = 0. 

V. Bayesian Blanking Nonlinearity Estimator (BNE) 
Fig. |2(a) and Fig. |2(b)| suggest that for highly impulsive noise behaviors (e.g., A = 0.001), the OBE 



shape resembles the BN shown in Fig. |4l where the received signal is simply blanked to zero when its 
absolute magnitude overpasses the threshold a. 




Fig. 4: . Blanking nonlinearity estimator 



Similarly to the SLE, the goal is to derive the optimum threshold foi" ^^N, in the MMSE 

sense. In this case the estimation error is 

sbn = hBNix, n;a) = X- gsNix + n; q), (19) 
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where hBN{x,n;a) is expressed by 

{— n , |x + n| < a 
(20) 
,\x + n\ > a 

By derivations similar to those of the SLE, detailed in Appendix C, the optimum a^^ is given by the 
solution of the following fixed-point equation 

j-,{mse) , N m=0 \ x m/ 

" = ^BN W = — • (21) 

^ m! ((t2 _^^2 )i/2e 

Although the fixed point problem admits a unique solution, as detailed in Appendix C, i^BN i^) 

a contraction mapping [21J, as it can be observed by Fig. [13] Consequently, F-^^^'^\a) in (|2TI) is not an 

attraction for the iterative algorithm Al, and an iterative algorithm that converges to the fixed point is 

A2: Iterative algorithm for optimal BN threshold 
set ao > 0, < ;U < 1, n = 0; 
while \F^'^''''\an) - a„| > e and n < UMax 

an+i = an + fJ-{an - F^'^'"'\an)); 

n = n + 1; 

end 

(mse) 

set = "n- 



In algorithm A2, the starting point oq has to be greater than zero (for instance ao = 2cr^) to avoid 
the useless solution a = 0, e represents the accuracy that is requested for the solution within a maximum 
number of iterations nuax, while ^ controls the speed of convergence. 

VI. Maximum SNR (MSNR) Estimators 

Any estimator x{y) is in general obtained as a non-linear transformation of the observation y = x + n. 
Anyway, it can be always decomposed as a scaled version of the input, plus a distortion term Wy, as 
expressed by 

X = g{y) = kyy + Wy, (22) 

where 

ky = EY{xy}/EY{y^} (23) 
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is the linear regression coefficient that grants the distortion term is orthogonal to the input term, i.e., 
E{ywy} = 0. In |[20l . |[23l it is suggested to exploit a similar approach to define a MSNR estimator for 
a complex Gaussian process (an OFDM signal) corrupted by an impulsive Gaussian mixture noise. Note, 
however, that due to the presence of the impulsive noise n, the non-linearity input y = x -\- n does not 
contain only the useful information x. 

Thus, in order to define a meaningful SNR ll20l . it is more convenient to express the estimator output 

as 

x{y) = g{y) = kxX + Wx, (24) 

where 

kx = EYx{9{y)x}/Ex{x^} = ExN{g{x + n)x}/Ex{x^}, (25) 

is the partial linear regression coefficient that grants the distortion noise Wx is orthogonal to x, as 
summarized by Exw{xWx} = 0. Although in general kx 7^ ky, when the inputs are both zero-mean 
Gaussian and independent it holds true that kx = ky |[26ll . By this approach, proceeding as suggested in 
|[20l . the SNR is expressed by 

o.rp_ klEx{x'} J^lEx{x^} 1 .... 

Ew{w^} Ej.{x^}-k^xEx{x^} EY{g{yr}/k^,al-l' ^ ' 



where the second equality in (1261 ) is granted by the uncorrelation between the useful part and the distortion 
noise. Clearly, the optimum non-linear estimator, in the maximum SNR sense, is expressed by 

xsNRiy) = argmai^[SNR] = argmm EY{g{yf}/ExY{9{y)x} ■ (27) 

In the problem at hand, taking into account that n is distributed according to ([3]l, the computation of kx 
in the denominator of (l26l) and (l27l) . can be obtained by 

_ ExN{g{x + n)x] _ r, ExNAgjx + n)x] 

where the subscript Nm means that the expected value is computed with respect to the m-th Gaussian 
pdf G{n; u^J associated to the Class-A Gaussian mixture. This fact suggests that the constant kx can be 
expressed as the weighted sum of other constants 

, _ ExNAajx + nm)x} 

which can be interpreted as the gain associated to the (virtual) useful components at the output of the 
non-linear estimator, when it is separately driven by the (virtual) input ym = x + Um- Note that, each 
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virtual input i/m is the sum of two zero-mean independent Gaussian random variables, thus, it is also 
zero-mean Gaussian distributed with variance = aj^ + cr^- 
In this virtual set-up, (see Theorem 1 in |26|), it holds true that 

, _ ExNAajx + nm)x} _ _ EYAg{ym)ym} 



where the expectation on the right-hand side of (1301 ) involves a single-folded integral, which is much 
simpler to compute than the two-folded integral in its left-hand side. Moreover, last equality in (l30l) 
is attractive because it is possible to exploit widely known results for the output of several non-lineai^ 
devices (such as the BN and the SL) driven by Gaussian inputs ll22l . ll27l - |[30l . 

Similarly, it is straightforward to derive that the average estimator output power, in the numerator of 
(l27l) . can be expressed as 

oo 

= EY{g{yf} = ^mEYAaHym)}, (31) 

m=0 



where also right-hand side of 0T| ) can exploit results widely available in the said technical literature for 
non-linear devices driven by Gaussian inputs. 

However, despite the above simplifications, the solution of the functional optimization problem in (l27l) 
is not easy, and consequently the derivation of the optimum estimator in the MSNR sense, is still an open 
problem. Conversely, it is possible to exploit (l30l) and (|3T]) . if g{y) is constrained to belong to families 
of suboptimal estimators, 

xxx{y) = gxx{y]a), (32) 

where XX stands for either the SL or the BN, and a is a scalar parameter that univocally specifies 
gxx{-, In this case the problem reduces to a classical optimization with respect to the scalar parameter 
a, where the optimum thresholds in the MSNR sense, for either the SL or the BN, are expressed by 

(snr) 



arg mm 



i?y{4x(y)}/(e^')' ■ (33) 



Thus, taking into account that the logarithm is a monotonic function, i.e., it does not change the position 
of extreme values, the optimum aj^x ^^"^ t)e obtained by solving 

A/iogM^xMUo (34) 

which leads to 

' ^ii;H4x(y)}-T7lol:4^^ = o- (35) 



EY{xl^{y)]da k^^^da 



13 



When the non Unear device is the SL gshiv', a) of Fig. [3l it is straightforward to derive ll22ll that 



, (SL) _ Ey^igsLiym; a)ym} _ f[ a 



y,m 



and, consequently, 



m=0 



a 



(36) 



(37) 



where the error function is defined as erf(x) = (2/^/7r) Jq e dt. Analogously, when the non-linear 
estimator is the BN gBNiv, ct) of Fig. |4l by standard integration techniques it is possible to prove that 

/,(BN) _ EYA9BN{ym;a)ym} _ ,(SL) ^ a 



a, 



y,m 



TT a, 



y,m 



and, consequently. 



^,(BN) _ jlCSL) 



a 



Similarly, plugging gshiy, a) and gB^iy, a) in (l3TI ). by standard integration techniques for Gaussian pdfs, 
it is straightforward to derive 

a \ 12 a 



(38) 



(39) 



Eyixlt^iy)} = E /^"ic^^.r 



m=0 



erf 



TT a. 



y,m 



and 



Eriildy)} = m^Uy)} + E 



m=0 



1 - erf 



(40) 



(41) 



Note that (|37] ). ( [39l) . (l40l ) and (|4T]) are different from the corresponding equations in |[20l . ll23l . because 
this paper deals with SL and BN of real random variables, while II20I and 11231 consider the SL and the 
BN of the envelope of complex random variables. However, if the source is a complex Gaussian random 
variable, (|37] ). ( [39l ). (l40l) . and (|4T]) can be applied to SL and BN that separately operate on each Cartesian 
component of the complex source. 



Plugging ( [38] ) and ( |40l ) in (|35] ). after standard algebraic manipulation the optimal BN threshold a 



(snr) 
BN 



is the solution of the following equation 



/3„ 



m=0 



^{^iN(y)} 



l(BN) 



0. 



(42) 



Analogously, plugging (l37l ) and (|4TI ) in (l35l) . the optimal SL threshold a"^^^^ is the solution of the 

'2 V- Pr, 



following equation 



I -erf (-^)) 



m=0 



^{^iL(y)} 



"'X 



(SL) 



0. 



(43) 



14 



Note that while (l42l) has to be solved by root-finding numerical techniques |[2ll . (1431) could also be 
solved by iterative numerical approaches similar to Al and A2, by exploiting its equivalent fixed-point 
formulation 



« = nr'i^) = ^itf^ • — , — r- (44) 



It is interesting to observe that the fixed-point problem in (144] ) can be expressed as 

which highlights that the fixed-point MSNR solution is in general different from the MMSE solution 
in ([TSl l. However, when the optimal thresholds are sufficiently higher than the input standard deviation 
ay, the power of the distortion noise is quite low and £^{x|^} f« k'^'^\, and, by (|43] ). the two optimal 
thresholds are very close if /c^ « 1. This fact can be further appreciated, by observing that the optimal 
MMSE thresholds are obtained by 

axx =argmm Ee{e ] , (46) 
and that, according to (l24b . the estimation error can also be expressed as 

e = x-x{y) = {l-kx)x-Wx. (47) 
Thus, exploiting the orthogonality between x and Wx, an alternative expression for the MSE is 

EE{e^] = (1 - kxfal + EwAwl], (48) 
and, exploitingE'y{2;xx(y)} = [^f^^'']'^ a\ + Ew^iw^.}, it is possible to conclude that 

EEie"} = (1 - 2^^))^ + Ey{xl^{y)}. (49) 
This means that the MSE optimal threshold is obtained by 

= argmin[£;y{4x(y)}-24^^Vi] 



argmin [^y{x|x(y)}/4^^' - 2^x1 
arg mill [Sy{x|x(y)}/4^^'] 



(50) 



arg mm 



which is different from the MSNR criterion in (1331 ) due to the absence of the square-power in the 
denominator of the cost function. Thus, MMSE and MSNR approaches provide very close thresholds 
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when A;^^^' 1, which happens, for instance, when the cUpping threshold a is sufficiently higher than 
Oy, due to a\j ^ a\ (see also (O). 

Noteworthy, the alternative expression (|49l ) is very useful for the computation of the MSE, because it 
requests to compute only single-folded integrals, i.e., the estimator average output power ii^y{xxx(y)} 
by dlD (i.e., & and dH]), for BNE and SLE) and the gain kf^^ by ^ (i.e., ^ and for BNE 
and SLE). Conversely, according to ([T3l) or ( fT9l) . the MSE should be computed by 

oo 

Ee{(?} = ExN{h'^ix,n;a)} = ^ /3mExN,Ah^i^^n'm;a)}, (51) 

m=0 

which involves tedious double-folded integrals with respect to the signal and the noise pdfs. 

VII. Computer Simulations 

This section verifies by computer simulations the analytical results derived so far. The MSE and SNR 
performance of the SL and BN are computed for several sets of the Class-A canonical parameters A, 
T, afj and for several SNR values. All the simulated MSE and SNR are obtained by generating 10^ 
observed samples y in Q. The Middleton's Class-A noise has been generated by the toolbox 131]. The 
optimal thresholds for the MMSE SL and BN are obtained by Al and A2, using e = 0.01 and fi = 0.01. 
The thresholds for the MSNR SL and BN are obtained by Madab® numerical solutions of (|43]l and (l42l) . 
respectively. The series with infinite terms, which are induced in all the analytical results by the Class-A 
pdf in Q, have been approximated by considering only the first M = 50 terms (although, M G [10, 20] 
would be enough in most of the cases). 

It is possible to observe in Fig. [5] and Fig. [6l the dependence of the optimal SL and BN thresholds on 
the total SNR, which is defined as SNRtot = crj^/{(^t + o"/)- It is shown that the MMSE and MSNR 
(optimal) thresholds are similar for high values of SNRtot (i-c, when aj^ ^ o"^), where, consequently, 
the two criterion are almost equivalent. This is not the case, for low (and negative) values of the SNRtot, 
where the MMSE and MSNR thresholds tend to diverge. Moreover, it is to note that the MMSE and 
MSNR thresholds tend to be more different for the BN rather than for the SL. This fact is more evident 
when T = (Jt / (jj = I, i.e., when the noise power is equally split between the AWGN and the impulsive 
component. 

Fig. El-Fig. [12] let to better appreciate the sensitiveness of the SNR and MSE performance with respect 
to the thresholds values of both the SL and the BN, and their performance penalties with respect to 
the OBE. Note that in all the figures the minimum MSE and the maximum SNR, correspond exactly to 
the thresholds values predicted by the theory. Moreover, also the theoretical MSE and SNR, perfectly 
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match with simulation results. As already anticipated, in most of the scenarios the optimal MMSE and 
MSNR thresholds are almost equivalent, and consequently they provide almost the same MSE and SNR 
performance. This is not the case however, in highly critical scenarios, when the SNRtot is quite low 
or negative, and the AWGN noise power is not negligible, with respect to the impulsive noise power 
(either because T f« 1 or because the impulsive noise tends to a Gaussian noise when A ?a 1). All the 
figures also highUght that the OBE, as theoretically expected, always outperforms the SL and the BN in 
MSE. The SL and BN penalties however are not dramatic, as it was expected by the fact that the OBE 



shapes in Fig. |2(a)| highly resemble either the SL or the BN, for several values of the canonical Class-A 
parameters. Interestingly the OBE derived in this paper, although it is not the MSNR optimal estimator, 
it outperforms in SNR the MSNR (sub-) optimal BN and SL, in most of the cases. 

CONCLUSIONS 

This paper has derived the MMSE Bayesian estimator for a Gaussian source impaired by impulsive 
Middleton's Class-A interference. The estimator is directly extensible to any Gaussian-mixture noise. Two 
popular and sub-optimal estimators, namely the soft-limiter and the blanker, have been optimized both 
in a MSE and SNR sense, deriving also closed form expressions for their MSE and SNR. A theoretical 
link between MSE and SNR at the output of the estimator has been established, and scenarios when 
the MMSE and the maximum SNR criteria are (almost) equivalent, or different, have been clarified. 
The theoretical analysis and computer simulations have shown that at least one estimator, among the 
optimum soft-limiter or the optimum blanker, can be always used as a sub-optimum estimator with 
minimal performance loss with respect to the MMSE Bayesian estimator. The derivation of the optimal 
estimator in the maximum-SNR sense, as well as a closed-form expression for the MSE of the optimal 
Bayesian estimator, are still open problems for possible further research. 
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Appendix A 

OBE DERIVATION 

S{f) = ^ {^{x)} is used to indicate the Fourier transform (FT) of s{x), defined by j {s{x)} = 

+00 

/ s{x)e~^'^'"f^dx. It is also reminded that the FT of a Gaussian pdf is still a normalized Gaussian 

—00 

function, as expressed by 

_2 NT riff. ^2 \ rut.Ji. \lri((\.J2 



(52) 



where g\ 



'Xf ~ ^/i^'^^x)- Thus, by exploiting the convolution and derivative properties of the FT 
transform, the integral in ( fTOl) can be expressed in the frequency domain by 



f {Pxix) * fmix)} = 7 {Px{x)] 7 {fm{x)\ 

= J {xG{x-a\)] f {G{x;al,)] 
By observing that the first derivative of a Gaussian function is expressed by 



(53) 



df 



then (I53]) becomes 



J 



7 {Px{x) * fra{x)} 

Thus, by the duality property of the inverse FT 

Px{x)*Ux) = -^^T-^[f GU;a\j)G{f-alj)} 



(54) 



(55) 



1 ^ T-^[G{f;a\f)GU,al^f)} 
a\i:[G{x-al)*G{x-al)]. 



Exploiting the well known property that the convolution of zero-mean Gaussian pdfs is still a zero-mean 
Gaussian pdf, with a variance equal to the sum of the two single variances 



Px{x) * fmix) = -ax-^G{x;ax + cj^) 
Summarizing, equation ([TOl l can be expressed by 



'X 



-xG{x; ax + cr. 



(56) 



00 2/3 



(57) 



m=0 



which coincides with (llTI) . 
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Appendix B 
sle derivation 

By observing ([T4l ) it is clear that hsh{x, n; a) is continuous with respect to a G TZ^ and that 

X + a , X < —a — n 



Consequently, substituting (1581 ) in ([T6l ). the optimum value Osle obtained by the solution of the 



,\x + n\ < a 

a — X ,x> —n + a 

(mse) 



(58) 



following equation 

+00 —a—n 



+ 00 +00 



{x + a)fx{x)dxfN{'n)dn+ J J {a — x)fx{x)dxf]y{n)dn = 0, 

— oo — oo —ooa—n 

which, by standard integration of Gaussian density functions, is equivalent to 

+00 



(59) 



a 



NsLja) 
DsL{a) 



-(Q+n)V2cr 



fN{n)dn 



+ 00 



/ [l>(-^)+$(^)]/iv(n)dn 



(60) 



where ^{x) = 0.5 [l + erf (x/\/2 



By substituting in (l60l ) the expression of fN{n) in it is recognized, by means of ([T]), that the 
numerator 



oo Am 



m=0 



ml 



+ 00 +0O 

J G{n; a'^)G{a — n; ax)dn + J G{n;a'^)G{a + n;ax)dn 



, (61) 



contains two integrals that represent the convolution, and the correlation, of two zero-mean Gaussian 
pdfs. Due to the even symmetry of Gaussian functions, the correlation is equivalent to the convolution, 
and the result is another zero-mean Gaussian pdf, with a variance equal to the sum of the two single 
variances, as expressed by 



oo oo . . -ry^ /2(rT^ 4- rr^ 1 



m=0 



m=0 

oo 



= 2a E Nm{a). 

m=0 

As far as the denominator D^i^{a) in ( [60l ). it is possible to recognize that 



(62) 



oo Am 



+ 00 



m=0 



Ox \ Ox 



G{n; a'^)dn, 



which can be rewritten as 



-A oo 



m=0 



DsLia) = l + "—Y.^j [-^rf + erf (l^)] G(n; ai)dn. 
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Exploiting the odd and even symmetry properties of erf(-) and G(-;ct^), respectively, and by a change 
of integration variable, it is finally obtained 

Dsda) = 1 + g 4;;-^Te'^"^""'+"^'/'"'"erf(^)dz 

"-"^ +00 , (63) 

= 1 + V ^^^e-"'/2<7,=^" f e~P^ e-''^' erf (z)dz 

„ m! WTTcr„^ J ^ ' 

m,=U —00 

where p = \f2aoxl'^%i > ^'^'^ ^ = '^x/'^m- The last integral in ( [63] ) can be rearranged as follows 
+00 00 

J e-P'e~'"" erf {z)dz= J e'^^ e'^'^d {z) dz+ j e-^^ e-^^\ri {z) dz , 
—00 —00 

which, by inverting the integration variable in the first term, turns to be equal to 
+00 +00 

j e-^^'erf (z) [e'^^ - eP^]dz = -2 ^ e'^^'erf (z) s\nh{pz)dz. 
—00 —00 
Substituting (64) in (l63l) . and exploiting 2.8.10.2 in fil)\, vol. 2, pag. Ill], it is finally obtained 

= 1 - e-^ E 



m=0 

-A V- iillerff- 



In order to prove the existence and the uniqueness of the solution to the fixed point equation (fTSl) . note that 
the relative minima and maxima of F^^^^ (a), are associated to the zeros of F^"^ (a) = dF^^"^^ (a) / da, 



expressed by 



m=0 m=0 V2'^('^i+'^™) 



m=0 



(64) 



Thus, the relative minima and maxima, correspond to values a* that satisfy the following expression 

00 

a* = ^ fifi = Ft'\a*)^^^ > Fir\a*). (65) 

m=0 ^ ™ m=0 ^ " 

Taking in mind that F^'^'''\Qi) = iVsL(O) > and Fsl'"''(0) = 4a|7VsL(0)^ > 0, the last inequality in 
(refeq:65) means that all the relative maxima and minima of F^^^'^ (a) occur when F^^'^'^ (a) is below 
the angle bisector q, as shown in Fig. (TO] Thus, necessarily F^^'^\a) crosses the angle bisector at 
least once, before its first relative maximum. Actually, to prove the uniqueness of such an intersection 
and, consequently, the uniqueness of the fixed point solution, it would be sufficient (not necessary) to 
prove that F^^^^ (q) is a contraction mapping ||2T1 between the fixed point solution and its first relative 
maximum (e.g., |FsL"(a)| < 1 for aopt < a < a\), as illustrated for more clarity in Fig. [TOl Although 
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it is not easy to analytically prove this inequality, it has been verified by brute force computer tests, 
partially reported in Fig. [14] , that it is always satisfied by any parameter set {A, T, ax} and any SNR 
value, and, consequently, it is conjecture as true. 

Appendix C 
bne derivation 

By observing (l20b . it is possible to proceed with the same approach of Appendix B. Due to the fact 
that ( [20l) is a piecewise constant function with respect to a G TZ^, with a discontinuity in a = |x + n|, 
the Dirac's impulse function 6{x) can be exploited to handle the derivative in such a point of the BN 
characteristic. Thus, by differentiating with respect to a, it is obtained that 

= [h\n+x\+{x,n)-h\n_^^\-{x,n)]6{a-\n + x\) (66) 

= —{n + x)S{a — \n + x\), 
where ha±{x,n) = lim hBN{x,n;a) represents the limit from either the right (+) or the left (— ). 
Consequently, by direct substitution of (l66l ) 



ha{x,n;a)h^^'°^''{x,n;a) = —{n + x)h{x,n;a)6{a — \n + x\) 

= -(n + x) '"+"'^ 'I '"+"' ^ ' 6{a-\n + x\) 
= —{n + x) ^~^^2 "'^ — \n + x\) 
= ^" 2^ ^ 5{a — |n + a;|). 

By exploiting the equivalent expression of ([TT]) for the BN, the optimum q, in the MMSE sense, is 

obtained as the solution of 

+ CXD +00 

(n^ — x'^) fx (x) Jn {n')S{a — |x + n\)dxdn = 0, 

— oo — oo 

which can be rewritten as 

- - . . . .... 7 7<. - . . . .... . 0. 

~oo — oo — oo — n 

By exploiting the integral properties of the Dirac's delta functions, it is possible to recognize that 

+ 00 +00 

— (— n — a)'^]fx{—n' — a)fj\f{n)dn + J [n^ — {a — n)^]/x(a — n)f]\f{n)dn = 0, 

— oo 

which, by simple manipulations, and exploiting the even symmetry of fx{x), can be simplified to 

+00 +00 

a / fN{n)fx{n + a)dn = / nfN{n)[fx{a-n)-fx{a + n)]dn. (67) 
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By substituting Q in (|67]) . the solution of the previous equation can be expressed as the solution to the 
fixed-point problem 

oo +00 

E Pm J nfm{n)[fx{a-n) - fx{a + n)]dn 

m=0 —00 //TIN 

a = ; , (68) 

00 +00 ' ^ ^ 

J2 l3m J fm{n)fx{n + a)dn 

m=0 —00 

where the integral at the denominator is the convolution of two Gaussian, zero-mean pdfs, which generates 

another Gaussian pdf with a variance equal to the sum of the two variances. Consequently, 

00 +00 

J2 Pm J nfm{n)[fx{a-n) - fx{a + n)]dn 
" = 55 z > (69) 

where the integral at the numerator can be split as the difference of a convolution and a correlation 

integral. By defining for convenience the odd function Pm{n) = nfm(n), it can be observed that 

+00 

A'BN(a) = j nfm{n)[fx{a- n) - fx{a + n)\dn 

—00 



= gm{a) * fx{a) - Qmi-a) * fx{a) 
= '^9m{a) * fx{a) 
= 2j-i{G™(/)Fx(/)}. 
By the same arguments that lead to equation (l56l ) in Appendix A, (TTOI ) becomes. 



(70) 



and, consequently, the fixed-point equation is expressed by 

00 

E 

;i(mse)/ N _ m.=0 



" = ^BN («) = ^5 > (71) 

m=0 

which coincides with (|2T]) . 

In order to prove the existence of the solution to the fixed point equation (TtTI ). first note that the problem 
has a different behavior with respect to the SLE case: indeed, F^^{a) is monotonically increasing, as 

proved in the following. Lets for convenience express F-Q^{a) in (TtTI ). by obvious notation equivalence, 

00 

AT ( \ ^ "m'^ 

^BN («) = «G(a) = Q , . = a^55 , (72) 

m,=0 



and its first derivative by 



FjJJ,'°)(a) = G(a) + aG^^^°'\a), (73) 



22 



where G{a) > for any a. Thus, in order to prove F^^^'^\a) is monotonic, it is only requested to prove 
that G^^'°'\a) > for any a. To this end, lets express G^^'°'\a) as 

2 V ''I J 



oo oo 



(74) 



and note that the terms for m = / null out in the double series in (|74]). Thus, due to the symmetry when 
the index m is interchanged with I, equation (1741 ) can be rearranged as 

Mi,a)/ \ a ^ fbrnai-a-mbi biam-aibjn\ -^{-^+^) 

^"^^^^n,^.A — — + — k — r ■ ^'^^ 

Gv / m=Ol=m+l ^ i / 

By substituting in (1751 ) the value of a^, bm, and fe^ subsumed in (TtTI ) and (r72l ). it is readily derived that 

Or \ / m=0 l=m+l '^m r^i 

due to the fact that all the terms inside the double series are greater than zero. Thus, it is proved 
that i^BN*'^^(a) is monotonically increasing, and that it starts from the non acceptable fixed-point so- 
lution -Fbn*^''(0) = 0' additionally, by observing that cj^ = maj/A + crj 00, it is proved that 
asymptotically Fg^*^-* (a) > a, as verified by 



lim ^^i^= lim 1^° = lim 7^ f = Ywl = 2- 

m=0 

(76) 

It is also highlighted, without a detailed proof, that -Fg^ ('^) starts to grow slower than a from a = 0, 
as expressed by -Fgjj'"^(0) = G{0) < 1. Thus, (see Fig. [TS] ) it is guaranteed that a solution ( at least) 
exist for the fixed point problem. Moreover, with the support of extensive simulation trials (see Fig. [T6l ) , 
it is reasonable to conjecture, that the fixed point problem for the BN has always the structure shown in 
Fig. [13] This conjecture guarantees the uniqueness of the (non-trivial) solution to the fixed-pint problem. 
Note that, the trivial solution a = should be avoided as starting point of the iterative algorithm A2. 
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Fig. 5: Optimal SL and BN thresholds {A = 0.01, T = 0.1, aj. = 1) 
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Fig. 7: MSE curves for A = 0.01 and SNRtot = dB 
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Fig. 8: SNR curves for A = 0.01 and SNRtot = dB 



27 




2 4 6 8 10 

SL and BN thresholds (a) 



Fig. 9: MSE curves for A = 0.01 and T = 0.01 




Fig. 10: SNR curves for A = 0.01 and T 
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Fig. 1 1 : MSE curves for T = 0.001 and SNRtot = dB 
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Fig. 12: SNR curves for T = 0.001 and SNRtot = dB 
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Fig. 14: Fixed point problem fertile SLE, wiien A = 0.1 T = 0.01, = 1, for severai vaiues of tine totalSiVi?. 



